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' We review the predictions of the theory of Color Glass Condensate for gluon produc- 

, tion cross section in p(d)A collisions. We demonstrate that at moderate energies, when the 

fT^ ■ gluon production cross section can be calculated in the framework of McLerran-Venugopalan 

' model, it has only partonic level Cronin effect in it. At higher energies/rapidities corre- 

, sponding to smaller values of Bjorken x quantum evolution becomes important. The effect 

^ • of quantum evolution at higher energies/rapidities is to introduce suppression of high-pT glu- 

I ' ons slightly decreasing the Cronin enhancement. At still higher energies/rapidities quantum 
evolution leads to suppression of produced gluons at all values of pr ■ 

■ 

> '. I. INTRODUCTION 

X' 

^ • Recently there has been a surge of interest in particle production in proton-nucleus (pA) and deuteron- 

5^ I nucleus (dA) collisions at high energies. The interest was inspired by the new data produced by the dA 



program at RHIC [1-4], which should enable one to separate the contributions of the initial state effects [8] 
such as parton saturation [9-13] from the final state effects, such as jet quenching and energy loss in the 
quark-gluon plasma (QGP) [14-17], to the suppression of high transverse momentum particles observed in 
Au - Au colhsions at RHIC [5-7]. 

Saturation physics has been largely successful in describing hadron multiplicities in Au-Au collisions 
at RHIC [18]. It can also have important implications for the transverse momentum distributions [19] 
and particle correlations and azimuthal anisotropics [20]. It has been demonstrated [21] that saturation 
provides very favorable initial conditions for the thermalization of parton modes with the transverse momenta 
kx ~ Qs, where Qs is the saturation scale. The thermalization was also found [21] to approximately 
preserve the centrality dependence of total hadron multiplicities determined by the initial conditions [18]. 
Recent lattice results [22,23] show that the initial average transverse momentum (kx) of the produced 
partons is (fcx) ~ Qs, which makes the "soft thermalization" scenario preserving the initial centrality and 
rapidity distributions quite likely. Final state interactions, however, will undoubtedly modify the transverse 
momentum distributions at fcy < (1 3) Qs without introducing new momentum scale [21]. If the produced 
medium lives long enough, then high kx jets will be suppressed as well because of the jet quenching and 
energy loss [14-17]. 

The first d — Au data from RHIC show Cronin enhancement extending up to fc^ ~ 6 GeV around y ~ 
[1,3] whereas at slightly forward rapidity around ?/ ~ 1 no significant enhancement is seen [2]. The absence 
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of suppression indicates that final state interactions are indeed responsible for the effect observed in Au — Au 
collisions [5-7] in the same kr range. However, the non-universality of the ratios for the charged hadron and 
TT*^ spectra [1] indicate deviations from the independent jet fragmentation up to kx ~ 5 GeV. Similar non- 
universality in the same kr range was observed for A and K production [24], and in p,p and pion production 
[25] in Au—Au collisions. It remains to be checked if there is a statistically significant suppression of high kr 
charged hadron and tt'^ yields above the Cronin enhancement region (fcy > 6 GeV), and if this suppression 
depends on centrality of d — Au collisions. This question is of crucial importance for the interpretation of 
the spectacular effect observed in Au — Au collisions [5-7] because this is the kinematical region in which the 
independent jet fragmentation picture, and thus the perturbative jet quenching description, begin to apply. 

The first dA data from RHIC [1-3] give the ratio of the number of particles produced in a dA collision 
over the number of particles produced in a pp collision scaled by the number of collisions 

R'^{kT,y) = . (1) 



To understand the new data on R and what it implies for our understanding of high energy nuclear wave 
functions we are going to study here the expectations for i?*^"^ from saturation physics. Our approach will 
be somewhat academic: in this paper, we will not include explicitly all of the effects related to the fact 
that high-fcr of produced particles corresponds to rather large Bjorken x in the actual RHIC experiments 
at central rapidity - the effective Bjorken x of high-fcr (fcr > 5 GeV) particles observed at mid-rapidity at 
RHIC at y/s = 200 GeV is about a; ~ 0.1 which may be too large for the small-x treatment that we present 
here (see [26], but see also [27]). These finite-energy effects have to be accurately accounted for before we 
can compare our calculations to the data. Nevertheless, we feel that a better understanding of the qualitative 
features of hadron production within the saturation framework is a necessary pre-requisite for a complete 
theoretical description of high energy p{d) — A collisions. 

We assume that collisions take place at very high energy such that the effective Bjorken x is sufficiently 
small for all kx of interest. For simplicity we will analyze proton-nucleus collisions assuming that the main 
qualitative conclusions would be applicable to dA. Since we can not calculate Ncoii in a model-independent 
way, we will be using Eq. (32) for our definition of i?^'^, which is identical to Eq. (1) applied to pA collisions 
with a proper definition of Ncoii (see [28] for a discussion of uncertainties involved in theoretical evaluations 
of this quantity). 

The problem of gluon production in pA collisions has been solved in the framework of McLerran- 
Venugopalan model [12] in [29] (see also [30-33]). The resulting cross section includes the effects of all 
multiple rescatterings of the produced gluon and the proton in the target nucleus [29]. At higher energy 
quantum evolution becomes important [34-41]. In the large Nc limit the small- a; evolution equation can be 
written in a non-linear integro-differential form [35-38] shown here in Eq. (49). The inclusion of non-linear 
evolution [35 38] in the quasi-classical gluon production cross section of [29] has been done in [42] (see also 
[43,44]). The study of the resulting gluon spectrum and corresponding gluonic R^^ is the goal of this paper. 

The paper is organized as follows. In Sect. II A we discuss two main definitions of unintegrated gluon 
distribution functions: the standard definition (2) and the one inspired by non-Abelian Weizsacker- Williams 
field of a nucleus in McLerran-Venugopalan model (6) [12,13]. We argue, following [42,45], that the Eq. (6) 
is the correct definition of the unintegrated gluon distribution counting the number of gluon quanta. We 
proceed by analyzing fcx-dependence of the distribution functions. In Sect. IIB we prove the sum rules for 
both distribution functions given in Eqs. (12) and (13), which are valid in the quasi-classical approximation 
only. In the framework of McLerran-Venugopalan model [12,46] the sum rules insure that the presence of 
shadowing in nuclear gluon distribution functions in the saturation regime {kr < Qso) requires enhancement 
of gluons at higher kr {kr > Qso) reminiscent of anti-shadowing. This conclusion is quantified in Sect. IIC 
[see Fig. 3 and Eqs. (28) and (29)] and the differences between distribution fimctions are clarified [see Eqs. 
(23) and (24)]. However, as we demonstrate in Sect. IIB, the sum rules break down once quantum evolution 
with energy [35-38] is included. They turn into inequalities (19) and (20). This indicates that, while multiple 
rescatterings in McLerran-Venugopalan model only redistribute gluons in transverse momentum phase space 
conserving the total number of gluons in the nucleus [29], quantum evolution of [35-38] actually reduces the 
number of gluons in the nuclear wave functions at a given value of Bjorken x. 

In Sect. Ill we study the gluon production cross section in pA in the quasi-classical approximation [29]. 
In Sect. IIIA we show that the gluon production cross section calculated in [29] in McLerran-Venugopalan 
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multiple rescattering model exhibits only Cronin-like enhancement [47], as shown in Fig. 4 and in Eq. (37) 
(cf. [48,49]). In the corresponding moderately high energy regime the height and fcr-position of the Cronin 
peak arc increasing functions of ccntrality as can be seen from Eq. (38). In Sect. IIIB following [42] we 
point out that, surprisingly, the gluon production cross section in pA can be written in a fcr-factorized form 
(43) [9,50] with the unintegrated distribution functions defined by Eq. (2), physical meaning of which is 
less clear than that of Wcizsackcr- Williams ones (6). In Sect. IIIB wc also prove a sum rule (46) which 
insures that suppression of produced gluons at low kx (fer^Qso) demands Cronin-like enhancement at high 
kr {kr ^ Qso) in McLerran-Venugopalan model. The relative amounts of suppression and enhancement are 
different from the quasi-classical gluon distribution case of Sect. II. 

Multiple rescatterings of partons inside the nucleus are believed to be the cause of Cronin effect. Phe- 
nomenologically these multiple rescatterings are usually modeled by introducing transverse momentum 
broadening in the nuclear structure functions [51-54]. In Sect. Ill we demonstrate how an explicit pQCD 
calculation of these multiple rescatterings done in [29] yields us Cronin effect (cf. [48,49]). 

Sect. IV is devoted to studying the effects of nonlinear evolution (49) on the gluon production cross section 
in pA. In Sect. IVA we use the analogy to the case of gluon production in deep inelastic scattering (DIS) 
solved in [42,44] to include the effects of evolution (49) in the gluon production cross section in pA. The result 
is given by Eq. (53). By expanding the all-twist formula (53) we then study the effect of nonlinear evolution 
on the gluon production at the leading twist (Sect. IVB) and next-to-leading twist (Sect. IVC) level. In 
Sect. IVB we start by deriving gluon production cross section at the leading twist level (62). We then 
estimate the cross section for high kx {kx > fcgoom ^ Qs{y)) in the double logarithmic approximation (66) 
[50] and demonstrate that the corresponding Rp^ is approaching 1 at high kr from below, i.e., that BP^ < 1 
at kx » Qs{y) [see Eq. (75)]. We proceed by evaluating Eq. (62) in the extended geometric scaling region 
{Qs{y) < kx < kgeom) [45,55,56]. The resulting leading twist gluon production cross section (86) leads to 
further suppression of gluon production due to the change in gluon anomalous dimension [8] as shown in Eqs. 
(89) and (91). At very high energies when the gluon production in pp is also in the extended geometric scaling 
region {kx < fc|eom) the ratio R^"^ saturates at i?^^ ~ A~^^^, as follows from Eq. (92). The next-to-leading 
twist contribution to the gluon production cross section in pA is evaluated in Sect. IVC with the result given 
by Eqs. (98) and (99). One can see that the subleading twist term contributes towards enhancement of gluon 
production at high kx- However, in the kx region where the next-to- leading twist contribution dominates 
over higher twists it is small compared to the leading twist term. Therefore the positive sign of the higher 
twist term can not alter our conclusion of high-fcy suppression we derived by analyzing leading twist. To 
understand how all twists add up we study what happens to Cronin peak {kx ~ Qs{y)) at high energy in 
Sect. IVD. We find that the height of the Cronin maximum decreases with energy and eventually Cronin 
peak flattens out at the same level as the rest of R^^ at higher kx, which is shown in Eq. (113). In Sect. 
IVE we observe that inclusion of evolution only strengthens the suppression of R^"^ at low kx {kx <C Qs{y)) 
[sec Eq. (116)] which was observed before in Sect. IIIB in the quasi-classical case. In Sect. IVF we construct 
a toy model illustrating the conclusion of Sect. IV that quantum evolution [35-38] introduces suppression 
of iJP"^ at all values of kx (see Fig. 8). We demonstrate that quantum evolution not only suppresses R^"^ 
making it less than 1, but also turns R^^ into a decreasing function of collision centrality contrary to the 
quasi-classical expectations. 

We conclude in Sect. V by summarizing our results. 

II. A TALE OF TWO GLUON DISTRIBUTION FUNCTIONS 

A. Definitions 

There are two different ways to define unintegrated gluon distribution function of a proton or nucleus. The 
most conventional way relates it to the qq dipole cross section on the target nucleus via two gluon exchange. 
Here we arc going to use a similar definition relating the unintegrated gluon distribution to the dipole cross 
section on the nucleus (see Fig. 1). 
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FIG. 1. "Conventional" definition of unintegrated gluon distribution relating it to the gluon dipole cross section. 
The exchanged gluon lines can connect to either gluon in the dipole. 



The corresponding gluon distribution is given by (cf. [41,43]) 
2^ Cf 



as (27r) 



(2) 



where Nq (r , 6, j/ = In 1 /x) is the forward amplitude of a gluon dipole of transverse size r at impact parameter 
b and rapidity y scattering on a nucleus [35,42]. We denote by k the transverse components of the four- vector 
k, and by kx its length. The definition of Eq. (2) is inspired by fcr-factorization and is valid as long as one 
can neglect multiple rescatterings of the dipole in the nucleus. By using Eq. (2) in the saturation region 
where higher twists (multiple rescatterings) become important one implicitly assumes that there exists a 
certain gauge in which the qq dipole cross section on a nucleus is given by a two gluon exchange interaction 
between the dipole and the nucleus and the interaction shown in Fig. 1 is literally all one needs to obtain 
the correct dipole cross section. It is not clear at present whether this is the case and such a gauge exists. 
Therefore the gluon distribution given by Eq. (2) does not give one the number of gluons in the nuclear wave 
function in the saturation region. The applications of the definition (2) will be clarified later. 

Another definition of unintegrated gluon distribution literally counts the number of gluons in the nu- 
clear wave function. To construct it in the quasi-classical limit of high energy QCD given by McLerran- 
Venugopalan model [12] one has to first find the classical gluonic field of the nucleus in the light cone gauge 
of the ultrarelativistic nucleus (non-Abelian Weizsacker- Williams field) and then calculate the correlator of 
two of such fields to get unintegrated gluon distribution function (see Fig. 2). 




FIG. 2. Definition of unintegrated gluon distribution in McLerran-Venugopalan model. 

The non-Abelian Weizsacker- Williams field of a nucleus has been found in [13], leading to the following 
expression for the corresponding gluon distribution [13,29] 



^'^'^ix,kf) 



2 7r2 



4Cf 



as (27r)3 



d^d^e-'!^-^ Tr(A^'^(0) • A^'^(r))) 



where 



(3) 



(4) 



with p the atomic number density in the nucleus with atomic number A, T{b) the nuclear profile function 
and A some infrared cutoff. 
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Generalizing Eq. (3) to include non-linear small-x evolution in it [35-38] is rather difficult. However, the 
problem of including small-x evolution has been solved for structure function and for the gluon production 

cross section in DIS [35,42]. Inspired by those examples wc conjecture that replacing the Glauber- Mueller 
[46] forward gluon dipole amplitude on the nucleus by its' fully evolved expression to be found from the 
nonlinear evolution equation [35,38] 

l_e-2:''3?oin(i/rTA)/4 ^ Ncir^h^y) (5) 

would give us the unintegrated gluon distribution function of a nucleus in the general case: 

<f>'^'^ix,e) = J dHcfre-'^-^ 1 NG{r,b,y = Inl/x). (6) 

Similar expression for gluon distribution was obtained earlier in [45] . 

An important observation concerning the two gluon distributions presented above has been made in [42,43]. 
It was shown that, while the Weizsacker- Williams gluon distribution of Eq. (6) indeed has a clear physical 
meaning of counting the number of gluons [13], it is the gluon distribution inspired by fcT-factorization and 
given by Eq. (2) that enters gluon production cross section in pA collisions and in DIS [42,43]. More precisely, 
the gluon production cross section including the effects of multiple rescatterings and quantum evolution in 
it can be reduced to a fcr-factorized form [9] with the unintegrated gluon distribution of a nucleus given by 
Eq. (2) [42]. The authors can not offer any simple physical explanation of this paradox. Nevertheless we 
keep both distributions in the discussion below keeping in mind that the first one is more relevant to particle 
production in pA. 



B. fcr-dependence: General Arguments 

Both definitions of unintegrated gluon distribution (2) and (6) have the same high-ZcT asymptotics in the 
quasi-classical approximation (see e.g. Eq. (3)), which reads 

Mx,k^) = 0A = A(t>N{x,k^) = A 1, kr ^ 00, (7) 

TT k 

where the index A (N) denotes gluon distribution in a nucleus (nuclcon). Therefore the distributions are 
equivalent at the level of leading twist, i.e., as long as we include only a single rescattering in the dipole 
amplitude Nq- 

In the quasi-classical case of McLcrran-Vcnugopalan model both gluon distributions obey a sum rule which 
we are going to prove here for <j). From Eq. (2) one can easily infer that 



/ d^kcj>A{x,e) = / dH [Vl NG{r,b,y = \nl/x)) 

J as (27r) J 



(8) 



r=0 



At very small rx the dipole cross section Nq in McLerran-Venugopalan model goes to zero as r|, (color 
transparency [57]) with the coefficient in the front proportional to A^/^. One can see that this is explicitly 
true for the Glauber-Mueller expression for the dipole cross section Nq [46] 

NG{r, b,y = 0) = 1- e'^^'^'o in(i/rTA)/4_ 

For Nq from Eq. (9) we observe that 

lim (V^ NG{r,b,y = 0) - A^/^V^ nG{r,b,y = 0)) = (10) 

where nc is the gluon dipole cross section on a single nuclcon obtained from Eq. (9) by expanding it to the 
lowest non-trivial order and putting ^ = 1. Remembering that 

/ dH = ^2/3 f d% (11) 
J A Jn 
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we conclude from Eqs. (8) and (10) that in the quasi- classical approximation (see also [29]) 

J (fkMy = 0,k^) = A J d''kcf>N{y = Q,k^). (12) 

Similarly one can show that the Weizsacker- Williams gluon distribution in Eq. (6) obeys the same sum rule 
in the quasi-classical approximation 

j d^kcj^'X'^iy = 0,F) =Aj d'k^^'^iy = 0,k'). (13) 

However, the sum rules of Eqs. (12) and (13) break down when the non-linear evolution with energy 

[35,38] is included. To sec this wc first note that for very small rr one can use the expression for No given 
by the double logarithmic approximation [50,56,45] (sec Sect. IV for details on similar calculations) 

i^GKTT y>,o,y) (2a,y)3/4 ^^^> 



with 



as = . (15) 



Similarly for the proton amplitude no we write 



2 a2 1 1/4 1 

where now the scale characterizing the proton is given by 

= 4nal^ (17) 

Op 

with Sp the cross sectional transverse area of the proton. Employing the fact that QIq = A^^^A"^ we can 
easily see that the amplitudes in Eqs. (14) and (16) do not satisfy the condition of Eq. (10) invalidating the 
sum rule. In fact using Eqs. (14) and (16) in Eq. (10) gives an inequality 

lim (vl NGir,b,y = Inl/x) - A^/^Vl nG(r,b,y = Inl/x)) < 0. (18) 
Eq. (18), together with a similar equation for Nc/rf' turn the sum rules of Eqs. (12) and (13) into inequalities 

j d'^k(l)A{x,k^) < A j d'^k(f)N{x,^) (19) 

and 

j d^k4>J'^{x,k^) < A j d^k^^'^{x,k^), (20) 

where the equality is achieved only in the quasi-classical limit. We conclude that while multiple rescatterings 
of gluons in McLcrran-Vcnugopalan model preserve the total number of gluons in a nuclear wave function at 
a given rapidity y, the quantum evolution tends reduce the number of gluons in the wave function via gluon 
mergers [9]. 

To study nuclear modification of the gluonic wave functions let us define the unintegrated gluon distribu- 
tions ratios as 

R^{x,k^) = 'f'^^'^'^'l and <^(a;,fc^) = '^^.li'^'-l . (21) 
'- ^ A^N{x,e) ^ ^ '- ^ A0]^>^(x,fc2) ^ ^ 
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The sum rules of Eqs. (12) and (13) imply that, in the quasi-classical approximation, if at some fcr the 
distribution function (j)^^ {y = 0, fc ) is smaller than A(j)^^{y = 0,k^), then at some other kr it should be 

bigger than A(p^^{y = 0, fc^). Using the definitions (21) one concludes from Eq. (13) that iiR^^{y 0, fc^) 
is below 1 at some kr it is bound to go above 1 at some other kx (for the same value of x). From Eq. (7) 
we can conclude that 

RA{y = 0,k^), Rj'^iy = 0,k^) ^1, kr ^ oo. (22) 

At the same time, when kr «C Qso the saturation effects become important driving (j)^^ [y = 0,fc^) below 
A(j)^^{y = 0,fc^), or, equivalently, making R^^{y = 0,k^) < 1. Therefore, due to the sum rule of Eqs. 
(12) and (13), somewhere in the region of kT>QsO the ratio Rj^iy = 0,fc^) should go above one. which 
corresponds to enhancement or broadening of the fcy distribution of gluons inside the nucleus. The same 
broadening argument applies to (j)A{y = 0,fc^). We have therefore proved that for both gluon distribution 
functions calculated in McLerran-Venugopalan model the effects of saturation and the sum rule (12), (13), 
while making RA{y = 0, fc^) < 1 in the infrared, also require an existence of a A:T-region where RA{y = 0, fc ) 
is above 1. This conclusion will be quantified in the next Subsection. 

The above argument docs not apply to the shadowing ratios i?^(.x, fc^) and R^^ {x, k^) when the effects 
of quantum evolution are included. The sum rules (12) and (13) are replaced by inequalities (19) and (20) 
which only require a reduction of the overall number of gluons in the nuclear wave function at a given rapidity 

y- 



C. fcT-dependence: Quasi-Classical Approximation 

To investigate the A:T-dependence of the unintegrated nuclear gluon distributions (j)^^ {x,k^) and (j>A {x, fc^) 
more quantitatively and demonstrate the difi^erences of the two distributions let us study them in McLerran- 
Venugopalan model [12,13]. For that we take the gluon dipole amplitude in the Glauber-Mueller approxi- 
mation [46] of Eq. (9). The high-fc^ asymptotic for both (j)^^ {x,l^) and (j)A{x,l^) is given by Eq. (7). 

Inside the saturation region (fcy Qso) one has 



<t>{x,kf) 



—2-, kT<^QsO, 



as (27r)2 g2„' 



(23) 



and 



6^^(a;,fc2) - ^^^^^ 



as (27r)2 kr ' 



In— ^, kr^QsO, 



(24) 



where we assumed for simplicity that the nucleus is cylindrical in which case its cross sectional area is 
Sa = -k'R? and Qso is given by Eq. (4) with pT{V) = A/Sa- 



Qlo = 



47r al A 



cylindrical nucleus. 



(25) 



In Eqs. (23) and (24) the difference between the two gluon distribution functions becomes manifest: 
4>^^{x,k^) keeps increasing (though only logarithmically) as kx decreases, while 4>Aix,k^) turns over and 
goes to zero in the infrared. Still for both distribution functions the ratio Ra goes to zero as kr ^ since 
to obtain Ra one has to divide Eqs. (23) and (24) by A0jv(x,fc^) from Eq. (7). The sum rules (13) and (12) 
require a region of enhancement {Ra > 1) at kx^ Qsa- To sec that the enhancement really happens one has 
to calculate the next-to-leading twist correction to the high-^T asymptotic of Eq. (7). This technique has 
been applied previously for quark production in [58]. One obtains 



i(x,fc2) 



C'f Sa QIq 
as (27r)2fc2 



1 + 2 



Q: 



3 +, 



ki 



oo, 



(26) 



and 



7 



as (27r)2fc2 

with 7 the Euler constant. For the ratios Ra's this imphes 



fcr oo, 



+ 27-3 + 



and 



r,WW 

St A 



(27) 



(28) 



(29) 



Therefore the ratios of gluon distributions approach 1 from above for both distribution functions at large 
kr- This of course indicates the presence of high-fcy enhancement. 

Ra 

2 . 5 

2 

1 . 5 

1 

. 5 



k/Qs 

FIG. 3. The ratio Ra of unintegrated gluon distributions in the nucleus and in the nucleon. The thin line represents 
the Wcizsackcr- Williams gluon distribution [Eq. (6)] while the thick line correspond to the more conventional one 

inspired by fcr-factorization [Eq. (2)]. 




Qualitative plots of ratio Ra for both distribution functions in McLcrran-Vcnugopalan model arc shown 
in Fig. 3. The thin line corresponds to Wcizsackcr- Williams gluon distribution (p^^ (x, k^) while the thick 
one represents the fc^-factorization distribution (j)A{x,k^). One can see that in accordance with Eqs. (23) 
and (24) the distribution 4>A{x,k^) goes to zero faster than (j)^^ {x, k^) as fc^ ^ in Fig. 3. In agreement 
with Eqs. (29) and (28) Ra for the distribution (pA{x,k^) has a stronger high-Zcr enhancement than R^^ 
for the distribution ^^^(a;,fc^). 

Finally, let us point out that the function Ra (R^^) shown in Fig. 3 will be modified when quantum 
evolution is included. Due to the inequalities of Eqs. (19) and (20) the total number of gluons will decrease. 
As we will see below in Sect. IV the effects of quantum evolution is to introduce suppression of gluons at all 
kr- 



III. QUASI-CLASSICAL APPROXIMATION: CRONIN EFFECT ONLY 



A. Gluon Production in pA 



The problem of gluon production in proton-nucleus collisions in the quasi-classical approximation 
(McLerran-Venugopalan model) has been solved in [29] (see also [30,33,31,32]). For a quark-nucleus scatter- 
ing the production cross section reads [29] 
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dgP^ 
dPk dy 



JdHd^x 



d'y 



1 asCp x-y 
(27r)2 7r2 iV 



which then has to be convoluted with the hght cone wave function of a quark in a proton. The saturation 
scale Qlo in Eq. (30) is given by Eq. (4). As was shown in [29] in the approximation when the logarithmic 
dependence of exponential factors in Eq. (30) on the transverse size is neglected, ln(l/a;TA) w x^, the xi_ 
and 2/_L integrations in Eq. (30) can be done exactly yielding 



dgP^ 
d^k dy 



as Cf 



d^b 



1 2 _j.2 ,^2 1 
1 p / VsO -I 



k' k' 



In- 



Ei 



(31) 



where Ei(x) is the exponential integral. Our goal is to construct the ratio of the number of gluons produced 
in a pA collision over the number of gluons produced in a pp collision scaled by the number of collisions 



R''^{k,y) 



<Pk dy 
<Pk dy 



(32) 



In the same approximation in which Eq. (31) is derived the gluon production cross section in pp scaled up 
by A is given by 



daPP 
d^k dy 



as Cf 



2i, QsO 



d^b 



(33) 



which can be obtained for instance by taking the kr/QsO ^ 1 limit of Eq. (31) and using the fact that 
QIq A^/^. For a cylindrical nucleus the impact parameter b integration would just give a factor of Sa- 
Using Eqs. (31) and (33) in Eq. (32) we then obtain 



RP^ikr) 



sO 



In 



Q 



Ei 



sO 



(34) 



The ratio RP^{kT) is plotted in Fig. 4 for A = 0.2 Qso- It clearly exhibits an enhancement at high-fcr typical 
of Cronin effect [47]. Similar conclusions regarding formula (30) have been reached earlier in [48]. 

It is worth noting that expanding RP^{k) from Eq. (34) in the powers of Qso/kr ("twists") yields a series 
with only positive terms 



RP^{kT) = 1 + 2^+6^ + 24^ + ... = f^n! 

— — — n=0 ^ — 



(35) 



The series (35) is divergent, but it is Borel resummable with the sum given by Eq. (34), though not all terms 
in Eq. (34) can be reconstructed by Borel resummation procedure. 

To establish whether inclusion of the correct transverse size dependence in the exponents of Eq. (30) 
would change the conclusion about Cronin effect let us study the high-AiT asymptotic of Eq. (30). A simple 
calculation yields 



dgP^ 
d?k dy 



Ois Cf 



d^^ 



In 



4A2 



+ 27- 1 + 



4k^ 



6 In 



-8(4-37) +29 + 247^-647 



kr — > 00. 



(36) 



For a cylindrical nucleus, keeping only the leading logarithmic (In^) terms in the parentheses of Eq. (36) 
we obtain 
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Rr'AikT) = l + l^ln^ + ..., kT^oo (37) 

indicating that R^^ approaches 1 from above at high kx, which is typical of Cronin enhancement. We 
therefore conclude that in the framework of the quasi-classical approximation employed in [29] the ratio R^"^ 
is less than 1 at small kx < QsO and has Cronin enhancement at high kr > QsO- 




k/Qso 

FIG. 4. The ratio R^^ for gluons plotted as a function of kr/QsO in the quasi-classical McLerran-Venugopalan 
model as found in [29]. The cutoff is A = 0.2 Q^. 



As can be seen from Eqs. (34) and (37), the position of Cronin maximum is determined by the saturation 
scale, such that /cmax = PQsO, where (3 is some weakly increasing function of \nQso/A. The height of the 
maximum is given by i?^'^(fcmax) — R^^{(3Qso)- Substituting kr = PQsO in Eq. (34) we observe that the 
height of Cronin maximum scales like 

R^^iP Qso) ~ In ^ + const ~ In^ -|- const'. (38) 

Since, for realistic off-central collisions A is replaced by the number of participants Npart, we conclude from 
Eq. (38) that in the quasi- classical approximation considered here the /cr-position and the height of the 
Cronin peak should increase with centrality of the pA collision. 



B. feT-Factorization 



Let us now show that it is possible to rewrite Eq. (30) in a /cT-factorized form [9,43,42]. Repeating the 
steps outlined in Sect. IV of [42] we first perform one of the transverse coordinate integrations in Eq. (30) 
rewriting it as 



1 asCp 



daP^ 



cPk dy 2n'^ tt 



dHd'^ze-'^-^ 



'■z-k . 1 
2 1 — In 



21.2 



z-^k 



ztA 



NG{z,b,0), 



where NG{z,b, 0) is given by Eq. (9). Using the fact that Ng{z = 0,b,0) = we write Eq. (39) as 

daP^ 



d^k dy 



^--^^JdHd^zNaU,mvl (e-^-Mn-i^). 



(39) 



(40) 



Let us denote the forward scattering amplitude of a gluon dipole of transverse size r on a single nucleon 
(proton) integrated over the impact parameter b' of the dipole measured with respect to the proton by 
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d^b' nG{r,y ,y = 0) = TTa^r^ In 



rr A' 



(41) 



Eq. (41) is obtained by expanding Eq. (9) at the leading order and taking A = 1. It corresponds to the 
two gliion exchange interaction between the dipole and the proton. In the quasi- classical Glauber-Mueller 
approximation in which Eq. (9) is derived each nucleon exchanges only two gluons with the dipole [46,13]. 
Therefore Eq. (41) is a natural reduction of Eq. (9) to a single nucleon case. 
With the help of Eq. (41) we rewrite Eq. (40) as [42] 



(Pk dy 



= ^-l3 4 / d'Bd%d'zylnG{z,b-B,0)e-'!^-^VlNG{z,b,0). 



(42) 



Now B is the impact parameter of the proton with respect to the center of the nucleus and b is the impact 
parameter of the gluon with respect to the center of the nucleus as shown in Fig. 5. 



proton 




gluon 



FIG. 5. Gluon production in pA collisions as seen in the transverse plane. To make the picture easier to read the 
gluon is placed far away from the proton which is highly unlikely to happen in real life. 

Eq. (42) is the expression for gluon production one would write in the /cr-factorization approach [43]. 
To see this explicitly let us rewrite Eq. (42) in terms of the unintegrated gluon distribution function from 
Eq. (2). One easily derives 



dgP^ 
d^k dy 



2a 

~C'f k 



^ J d^q(t)p{q)(l}A{k- q), 



(43) 



which is the same formula as obtained in fcT-factorization approach [9,50,43]. (j)p is defined as unintegrated 
gluon distribution of the proton given by Eq. (2) with nc instead of Nq on the right hand side. Eq. (43) 
demonstrates that the gluon production cross section in pA can be expressed in terms of the gluon distribution 
(2) in a rather straightforward way [42]. Somehow it is the distribution (2) and not the Weizsacker- Williams 
distribution (6) that enters Erj. (43). 

Eq. (43) demonstrates that, at least in the framework of McLerran-Venugopalan model, the multiple 
rescattering leading to Cronin enhancement in pA can be incorporated in the gluon distribution functions 
[29,33]. There is no clear distinction between the nuclear wave bmction effects and the Glauber-type rescat- 
terings in the nucleus. Anti-shadowing present in the gluon distribution function 4>a{M} as shown in Fig. 3 
simply translates into Cronin effect of Fig. 4 via Eq. (43). 

In the quasi-classical approximation of McLerran-Venugopalan model one can prove a sum rule for the 
gluon production cross section in pA similar to the sum rule we proved for gluon distributions in Sect. IIB. 
To prove the sum rule we note that Eq. (42) implies that 



/ 



d^kk 



daP^ 



d'^k dy 



as 27r2 



/ d^B d% [Vl nG{z,b- B, 0)] [V^ NG{z,b, 0)] 

J z=0 



(44) 



z=0 



For Glauber-Mueller Ng from Eq. (9) and for ug from Eq. (41) the following condition is satisfied 

^lim { [V^ UGiz, b-B, 0)] [V^ Ng{z, b, 0)] - [y^ nG{z, b-B, 0)] [V^ hgU, b, 0)] } = 0. (45) 
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The impact parameter integration in pA will give an extra factor of A^l^ as compared to pp. Together with 
Eq. (45) this gives 

in the quasi-classical approximation. 

Similar to the sum rule proved in Sect. II for gluon distribution functions, the sum rule (46) insures that 
if the quasi-classical gluon production cross section in pA collisions is, in some region of fcr, smaller than A 
times the gluon production cross section in pp than there should be some other region of in which their 
roles arc reversed. For defined in Eq. (32) that means that if, in some region of /c^, it is less than 1 
there must be some other region of A:t in which it is greater than 1. Of course the ^ factors in Eq. (46) 
make the quantitative amounts of suppression and enhancement very different from the ones dictated by, for 
instance, the sum rule of Sect. II. 

In the quasi-classical approximation for the gluon production in pA considered above EP^ is below 1 at 
kr'^QsO- Expanding Eq. (34) for fcy <C Qso we write 

7,2 

R^^ik) « < 1 if fcT « QsO. (47) 

Eq. (47), together with the sum rule (46) imply that there must exist a region of kx with a Cronin-like 
enhancement of gluon production, which is demonstrated by the full answer plotted in Fig. 4. 



IV. INCLUDING SMALL-X EVOLUTION: SUPPRESSION AT ALL Pt 



A. Including Small-a; Evolution 

As the energy of the collisions increases quantum evolution corrections become important. For produced 
particles with the same kx higher energy implies smaller effective Bjorken x meaning that the quantum 
corrections of the type a^lnl/.T should be resummed. These corrections can be resurnmed by the BFKL 
equation [34], which calculates the contribution of the hard (perturbative) pomeron. However, as energy 
increases multiple pomeron exchanges become important, resulting in a more complicated small-a; evolution 
[39,40]. In [35,38] an equation was constructed which resums multiple pomeron exchanges for a forward 
amplitude of a qq dipole scattering on a nucleus in the large limit. The forward amplitude N{r,b,Y) of 
a dipole of transverse size r scattering at impact parameter b and rapidity Y was normalized such that the 
total qqA cross section was given by 

al!^' = 2 J dHN{r,b,Y). (48) 
The evolution equation for N{r,b,Y) closes only in the large- A^c limit of QCD [38,39] and reads [35-37] 

Nix,„b, Y) = N{x^,A y = 0) '"(T)^ + ^ Z"" dy e"^ H^^Y^-y) 

f 1 1 1 

X / rf2a:2^-2i2-[2A^(^02,^+;^^i2>y)--^(^2>^+;^^i2>2/)-^fe2.^+;5^02>2/)]> (49) 

with the initial condition given by N{xQi,b, y = 0) taken to be of Mueller-Glauber form [46] in [35]: 

N{x^„bo,Y = 0) = 1 - e-2oiQSr''in(i/^oiTA)/4^ (50) 

where 
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In [42] it was shown how to resum the effects of nonlinear evolution of Eq. (49) for gluon production 
in DIS. In the quasi-classical approximation the gluon production in DIS is given by a formula similar to 
Eq. (30) [44] . That formula can also be recast in a kx factorized from of Eq. (42) [42] . As was proven in [42] 
in order to include quantum evolution (49) in Eq. (42) for DIS one has to make replacements. First, one 
has to replace NG{z_,b,0) in Eq. (42) by the forward quark dipole amplitude using the following expression 
valid in the large- A^c limit [42] 

NgU, b, y) = 2N{z, b, y) - N{z, b, yf, (52) 

where N{z_, b, y) is the forward scattering amplitude of a qq dipole on a nucleus evolved by nonlinear equation 
(49). Then one has to replace nf;(z, 6, 0) by na(z_,h,Y — y) evolved just by the linear part of Eq. (49) (the 
BFKL equation [34]). Here Y is the total rapidity interval between the projectile (virtual photon) and target 
nucleus in a DIS collision. The initial conditions for both Nq and no evolution are given by NG{z,b, 0) and 
nc {z, b, 0) correspondingly. 

Since both pA and DIS are considered here as scatterings of an unsaturated projectile (proton or qq pair) 
on a saturated target (nucleus) with the gluon production in the quasi-classical limit given by the same 
Eq. (42), we may conjecture that inclusion of quantum evolution (49) in gluon production cross section is 
done similarly for both processes. We therefore write 

daP^ Cf 



,2.. - j d^Bd%<fzVlnG{z,b-B,Y-y)e-'^-^VlNG{z,b,y), (53) 

a^A; dy as7r(27r)'^ J 

where Y is the total rapidity interval between the proton and the nucleus. Just like in DIS Ng in Eq. (53) is 

given by Eq. (52), where A'^ should be foimd from Eq. (49), while ug should be determined from the linear 
part of Eq. (53) (BFKL) with the initial conditions given by Eq. (41). Eq. (53) is exact if the proton is 
modeled as a diquark-quark pair [59] , in which case it would be identical to qq pair produced by a virtual 
photon in DIS. In general case Eq. (53) remains a well-motivated conjecture. 

Like in Sect. II the sum rule (46) breaks down once non-linear evolution [35,38] is included in the way 
shown in Eq. (53). Using the double logarithmic expressions (14) and (16) modifies Eq. (45) into 

hm { [Vl nG{z, b-B,Y-y)] [Vl Ng{z, b, y)] 

zt — ^■U 

- [V^ UG {z,b-B,Y- y)] [Vl TIG {z, 6, y)] } < (54) 
turning the sum rule of Eq. (46) into an inequality for the cross section from Eq. (53) 



ID ) into an inequality lor tne cross section 
/■ j2, ,2 daP^ ^ A f ,2 daPP 



(55) 



Again the effect of quantum evolution is to reduce the total number of gluons at a given rapidity, though now 
it is shown for the case of gluon production weighted by fcf.. Let us now study in detail how this suppression 
sets in for various regions of kr- 

In the following we are going to study effects of evolution equation (49) on the gluon spectrum and on EP-^. 
Our goal is to determine whether Rp^ preserves the shape shown in Fig. 4 with Cronin maximum and low-fcr 
suppression, or quantum evolution would modify this shape introducing extra suppression. Below we will first 
study the effects of quantum evolution at high-fc-r, ^t^Qs, showing that evolution introduces suppression 
[RP^ < 1) in that region. We will then proceed by studying the fate of the Cronin peak (fcr ^ Qs) as 
evolution sets in. We will show that the Cronin maximum will decrease with the onset of evolution and 
would eventually disappear. We will then argue that suppression persists for kr <C Qs when evolution effects 
are included. We will end the section by constructing a toy model summarizing our conclusions. 
To simplify the discussion we will consider a cylindrical nucleus for which Eq. (53) reduces to 



/ d'zVl 



d^kdy as7r(27r)3 ^ 



nG{z,Y-y) e-'^-2 Ng {z, y), (56) 



with Sp the cross sectional area of the proton. 
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B. Leading Twist Effects 



1. Leading Twist Gluon Production Cross Section 

We start by exploring the leading high-fc-r behavior of the gluon spectrum given by Eq. (56). At very 
high kT the integral in Eq. (56) is dominated by small values of zt- Therefore we can neglect the quadratic 
term in the evolution equation for Nq (49) leaving only the linear part - the BFKL evolution with initial 
conditions for a gluon dipole given by Eq. (9). The corresponding Feynman diagram is shown in Fig. 6. The 
solution of the BFKL equation is well-known and reads 

Naiiz^y)^ J ^^C^{zTQso)^e''''^^^^y (57) 

with 

x(A)=^(l)-i^(l-^)4v^Q), (58) 

as defined in Eq. (15) and Qso for a cylindrical nucleus given by Eq. (25). The coefficient is fixed from 
the initial conditions at y = given by Eq. (9). Then for small zt < 1/Qso 



^ ~ 2^ ^ 4" (n- m)! (2n- A)'"+i A 

ri=l m=0 ^ ' ' 

= y {^^'\{x + n, (2„_A)ln%V (59) 

Z^^ 4" n\ (A - 2n)"+i V ^ / V ^ J 

Similarly for the gluon dipole cross section on the proton we write 

naiz^y)^ j ^^ClizTkfe'^'^^^^y, (60) 

where the scale characterizing the proton A is given by Eq. (17). The coefficient is obtained by requiring 
that Eq. (60) reduces to Eq. (41) when y = 0. For zt < 1/A we derive 

(We have identified the non-perturbative scale characterizing the proton (17) with the infrared cutoff em- 
ployed earlier in Eq. (41).) 

D- proton 
BFKL 




BFKL 
nucleus 

FIG. 6. Gluon production in pA collisions at the leading twist level (see text). 
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Substituting Eqs. (57) and (60) into Eq. (56) and integrating over z_ yields [60,61] 



<Pk dy 



S„ 5.4 f dX d\' 2 2 ^ 2A+A' ( QsO Y f A 



4Q!^(27r)3y 27ri 27ri ^' r (2 - ^) V / V^t 

^ g2 a, [x(A) y+x(A') (Y-j/)] _ (-52) 



Eq. (62) gives the leading twist expression for the gluon production cross section in pA collisions and is 
illustrated in Fig. 6. 

The difference between Eq. (62) and Eq. (13) of [60] is in gamma-functions in the integrand. The difference 
manifests itself at the order of higher twists, where the gluon distributions used in [60] , if taken at y = and 
used in inverted Eq. (2) to obtain Up, would yield higher twist corrections (higher powers of r^) to the right 
hand side of Eq. (41), which should not be there in the two-gluon exchange approximation corresponding to 
y = limit. 



2. Double Logarithmic Approximation: Monojet Versus Dijet and First Signs of High-pT Suppression 

To derive the high-fcy behavior of the gluon production cross section in Eq. (62) we have to evaluate the 
integrals in it by saddle point method. When kr ^ QsO) A we approximately write 



X(A) 



Then the saddle points are given by 



^sp — 2 



2- A 



2as y 



ln(fcT/Qso) 



(63) 



(64) 



and 



Kp — 2 



^ 2as{Y-y) 
\n{kT/K) ' 



(65) 



Integrating over A and A' around the saddle points (64) and (65) in Eq. (62) yields gluon production cross 
section in double logarithmic approximation (DLA) [50] 



d?k dy 



DLA 



Cp Sp Sa Qso A^ 1 
as (27r)4 ^4 ^ 



y3 (r _ y)3 




xexpf2j2a«yln^ + 2-y'2a,(F-?/) In ^ 



(66) 



To understand Eq. (66) let us first construct gluon distribution function 



xGa{x, Q 



A2 



dk^ (l)^{x,k^) 



(67) 



in the same double logarithmic approximation [10]. Using Eq. (57) in Eqs. (2) and (67) we obtain in the 
double logarithmic approximation 



a,(27r)3 " (2S, ln(l/x))3/4 



lni/4 _Q_ ^ 

Qso ^2y/2a,\n(l/x)\n(Q/Q,o)_ 



(68) 
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Since the above gluon distribution is obtained in the DLA (large Q"^) limit of the BFKL equation, it can also 
be obtained by taking the small-a; limit of the DGLAP equation [62] . One can explicitly check that with the 
help of Eq. (68) and an analogous one for the proton gluon distribution xGp, Eq. (66) can be rewritten as 



dgP^ 
(Pk dy 



2 as 
Cpk^ 



xGp{x 



d 



xGAix = e-y,k^) + 



+ xGa{x = e-y,k^)-^xGp{x = e-^+^fc2) 



(69) 



Eq. (69) can be obtained directly by using Eq. (68) in Eq. (43) and assuming that the g-integration in 
Eq. (43) is dominated by the regions near q = and q = k [9,50]. As was shown in detail in [11], Eq. (69) 
can be reduced to 



d?k dy 



dx\ 

-Y + y X\ 



xiGp{xi,k^) xGa{x = e ^,k^) + 



+ 



— xGp{x = e-^+y,l^)xiGA{xi,^)\ -^pj—-^ 



(70) 



which is the standard dijet production cross section derived in coUinear factorization approximation (see e.g. 
[11]). (Of course one of the jets' momentum in Eq. (70) is integrated over.) Therefore, we have started with 
a single jet production cross section given by fcT^-factorized expression (62) with BFKL gluon distributions 
and demonstrated that in the large kx limit it reduces to the conventional dijet production cross section 
(70) given by collinear factorization with DGLAP-evolved structure functions.^ 

Before wc continue let us study RP^ given by the cross section of Eq. (66). The naive expectation for 
the high-fcT limit of the leading twist gluon production cross section would be that RP^ = 1. However, 
already at the level of approximation employed in Eq. (66) this is not quite the case. To see this let us first 
write down an expression for the gluon production cross section in pp collisions in the leading twist DLA 
approximation. It is obtained by replacing Qsq and Sa in Eq. (66) by A and Sp correspondingly. We obtain 



d^k dy 



C p Sp 



1 vv + ^/y^ 



DLA 



A4 

as (27r)4 fc4 ^^ymjY- 



y)3/4 



exp 



2^2a, ln^(y^+v/F^) 



(71) 



To calculate we note that since Sa = A^/^Sp one concludes from Eqs. (25) and (17) that QIq = A^/^K^. 
Using Eqs. (66) and (71) in Eq. (32) yields 



BP^{kT,y) 



1/4 



feT»Qs 




X exp 



2^2 a«y Win 



kr 



In 



kr 



where Qs 



Qso V A 

Qs{y) is the full energy dependent saturation scale, which reduces to Qso at y = 0. Defining 



(72) 



(73) 



'^We thank Al Mueller for encouraging one of the authors (Yu. K.) to verify this correspondence explicitly several 
years ago. 
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we rewrite Eq. (72) as 



exp 



L_ , Qso 



(74) 



where ^ < 1 for fcr > QsO, since Q% = A^/^A^ ^ A^. For large transverse momenta in question, kx Qs{y), 
the variable ^ is approaching 1 from below as is clear from Eq. (73). In the limit ^ — > 1 Eq. (74) becomes 



Vv - yy - y 



exp 



2;/2a,y(l-0 In ^ 



< (2 - exp 



2a,y(l-e) In^ 



exp 



2\/2a,y(l 



In^ 



< 1, fcr » Qs(y)- (75) 



We neglected 2 — ^ = 1 + (1— ^) in front of the exponent in Eq. (75) since the (1 — ^ correction to 1 in it is not 
enhanced by any parametrically large variables, such as y and InQ^o/A in the exponent. If RP^{^,y) from 
Eq. (75) is expanded in powers of (1 — this prefactor term would give subleading logarithmic corrections 
to the expansion of the exponent, which arc negligible in the DLA limit considered here. 

We conclude that RF^ from Eq. (72) is smaller than one. Since RF^{£,,y) in Eq. (75) is an increasing 
function of ^ and ^ is an increasing function of kr, we observe that RP^{kT,y) in Eq. (72) is an increasing 
function of fcy approaching 1 from below. This suppression is mainly due to the difference of the cutoffs 
in the logarithms of transverse momentum in the exponent of Eq. (72). The cutoff for the nucleus case is 
given by nuclear saturation scale, which is different from the appropriate scale in a single proton. The high 
momentum regions, where linear evolution equations work, are cut off from below by saturation scales, which 
arc different for different nuclei and for the proton. In this way, as we can see from Eqs. (72) and (75), 
saturation influences the physics at high kr as long as corresponding xsj is small. The effect of saturation 
is to introduce high-Zcr suppression. 

The suppression of Eq. (75) is a leading twist effect and is due to quantum evolution. In this sense it is 
similar to the suppression suggested in [8]. However, the suppression of [8] corresponds to a region of lower 
kx, where the double logarithmic approximation of Eq. (72) is not valid any more. There the suppression 
happens due to the change in anomalous dimension of the gluon distribution function, as we are going to 
discuss below. 



3. Onset of Anomalous Dimension: More High-px Suppression 

For the values of kr lower than considered above (but still much larger than Qso) the saddle point of 
A- integration in Eq. (62) shifts to a smaller value than given by Eq. (64). While in determining the saddle 
point of Eq. (64) we had to expand x(A) around A = 2, now we have to expand it around A = 1. There one 
writes 

X(A) « 21n2+Ic(3)(A-l)2 (76) 
obtaining the value of the saddle point 

As was suggested in [45] , the transition of the saddle point from the value given in Eq. (64) to the one given 
in Eq. (77) happens around 

fcgeom « Qs{y) %^ (78) 
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indicating the onset of geometric scaling regime [55]. Here in the double logarithmic approximation the 
saturation scale depends on energy as [56,45,63] 

Qs{y) ~ Qsoe^^-y. (79) 

The precise value of the scale fcgeom in Eq. (78) depends on the definition of the saturation scale and on the 
way one defines the transition between the double logarithmic and geometric scaling regions. For instance, 
if we define the transition by equating the saddle points of Eqs. (64) and (77) 

Aap = A:^, (80) 

we get at the point of closest approach (the two saddle point values are never equal to each other) 

When combined with the saturation scale from Eq. (79) this gives 

W «^ (82) 

which is slightly different from Eq. (78). At the same time, using the energy dependence of the saturation 
scale found in [64] in the fixed coupling case 

Qs{y) « aoe'-^^"^" (83) 

in Eq. (81) gives 

( \ \ 0.34 

-AV) 



Algeom « Qs{y) j ' (84) 

which is even lower than Eq. (82). A definition of the transition point difiierent from Eq. (80) would give 
slightly different estimates for fegoom. 

Nevertheless, the ambiguities in the scale /Cgeom notwithstanding, one can argue, as was done in [45], that 
there exists a large momentum scale fcgeom, which is parametrically larger than the saturation scale 

fcgeom > Qs{y)- (85) 

For kr > fcgeom there is no geometric scaling and the gluon production is well described by double logarithmic 
approximation described above resulting in RP-^ from Eq. (72). kr < fcgeom is the region of geometric scaling 
[45]. When kx^Qsiv) (saturation region) multiple pomeron exchanges become important leading to the 
saturation of structure functions [9] . For Qg (y) <kT< kgeom (extended geometric scaling region) multiple 
pomeron exchanges are not important yet and the gluon production cross section is described by the leading 
twist expression in Eq. (62) with the A-integral evaluated near the saddle point of Eq. (77) [8]. 

Performing the A and A' integrals in Eq. (62) in the saddle point approximation around the saddle points 
of Eqs. (77) and (65) correspondingly yields 



d?k dy 



CpSpSA QsoA' Ct In'/^^ 



LLA "«(27r)4 fc' ^/Tm^s{Y-yf/^{2cisYI^ 



X exp 



{ap-l)y + 2 



2a,{Y-y) In^ 



lu- 



14C(3)a.y 



(86) 



where 



ap — 1 = 2 as In 2 



(87) 
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is the BFKL pomeron intercept [34] and is well-approximated by the first term in the series of Eq. (59) 
for all physically reasonable values of A 



l + ln 



Q_ 
A 



sO 



l + -lnA 



(88) 



The superscript ^^"^ in Eq. (86) denotes the leading twist contribution. We assume that in the transverse 
momentum region where Eq. (86) is valid, Qs(t/) < fcy fcgeom, the gluon production cross section in pp 
collisions is still given by Eq. (71). This is a good approximation since if kT>Qsiy) ^ A the double 
logarithmic approximation of Eq. (71) should work. Using Eqs. (86) and (71) in Eq. (32) we obtain 



R''^{kT,y) 



2C^ 



1^1/4 



,,1/4 



Qso ^703) (2a.)i/4 ^+yF" 



X exp 



(ap - 1) y - 2 W 2 as y In — 



In^ ^ 



A 14C(3)«sy 



(89) 



To determine whether RP^{kT,y) in Eq. (89) is greater or less than 1 we first drop the slowly varying and 
constant prefactors in front of the exponent and write 



Qs(y)<kT<kg, 



Q 



sO 



exp 



(ap - 1) 2/ - 2 W 2 as ?/ In — - 



ln■ 



2 _kT_ 



A 14C(3)asy 



(90) 



keeping only parametrically important factors. To estimate the value of i?^'^ in Eq. (90) in the extended 
geometric scaling region Qs{y) ^kT< fcgeom we substitute kr = fcgeom into Eq. (90) with /cgeom from Eq. (78). 
The result yields an asymptotically small value 



RP^ikT,y) 



Qs(y)<kT<k, 



(91) 



'geom 



where we used A = 197 for gold nucleus. For other values of A and for other values of kr in the re- 
gion Qs{y)^kT<kgeoin One still gets exponential suppression for BP-^{kT,y)- Therefore we conclude that 
RF^{kT,y) < 1 in the extended geometric scaling region Q.s(?y) ^ ^ fcgcom- 

As can be checked explicitly, for sufficiently large nucleus (large A), RP^{kT, y) of Eq. (89) is an increasing 
function of k for (5« (y) < Ajt < fcgeom • As kx increases it should smoothly map onto BP'^{kT,y) of Eq. (72), 
which would approach 1 from below for asymptotically high kx- 

At very high energy the geometric scaling regions for the nucleus and the proton will overlap. Namely, 
the geometric scale for the proton /c|eom = k^eom/A^/^ will become larger than the saturation scale for the 
nucleus Qs{y) allowing for a region of kx where anomalous dimension sets in for gluon production both in pA 
and pp. ^ In this asymptotic region one has to estimate the A and A' integrals in Eq. (62) around the saddle 
point given by Eq. (77) (with A instead of Qso for the A' integral). Replacing Eq. (71) by the appropriate 
expression where the saddle points of A and A' integrals were given by Eq. (77) with A instead of Qso we 
obtain the following asymptotic expression at mid-rapidity (y = Y/2) 



BF\kT,y) 



-1/6 



Qsiv)<kT<kl,, 



exp 



ln2 *|f -In^-^ 
14C(3)asy 



(92) 



^The onset of anomalous dimension does not imply saturation and is still a leading twist effect. Therefore Eq. (56) 
in which no saturation in the proton's wave function was assumed is still valid in this region. 
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From Eq. (92) we conclude that in the extended geometric scahng region at asymptotically high energies, 
W"^ saturates to a parametrically small lower bound, R^^ ~ A~^/^, which is independent of energy and is 
a decreasing function of A, or, cquivalcntly, ccntrality. 

To conclude our discussion of high-/cT suppression at the leading twist level we note that, as was recently 
argued in [65], the running coupling effects in the BFKL evolution may modify the A-dependence of the 
saturation scale given by Eqs. (79) and (83), making Qs{y) almost independent of A at very high energy 
corresponding to large rapidity y. This would result in high-A;T suppression which would not disappear at 
any kr- That is EF'^{kT, y) would not approach 1 anymore at high kr- Instead one would have BF^{kT, y) ^ 

3 



C. Next-To-Leading Twist 

Above we have shown that the effect of quantum evolution (49) on the leading twist gluon production cross 
section in pA with kr > Qs{y) is to introduce strong suppression of Rp^. Here we would like to study the 
effect of evolution on the gluon production at the next-to-leading twist level. Below we are going to show that 
if one includes the evolution of Eq. (49) into the ncxt-to-lcading twist correction to Eq. (62) it would start 
contributing towards enhancement of R^^ at high-ZcT- This appears to indicate that multiple rescatterings 
always tend to enhance gluon production at high-fcr- As we will argue later the effect of quantum evolution 
is much stronger. It dominates at high energies leading to overall suppression of Rp^. 

A perturbative solution of Eq. (49) was constructed in [36] giving the forward amplitude of a qq dipole 
scattering on the nucleus as an expansion in powers of rrQsiy) 

N{r,b,y) = N,ir,b,y) + N2{r,b,y) + ..., (93) 

where the leading behavior of the nth term in the series is Nn{r,b,y) ~ {fTQ s{y))"' ■ To find the next-to- 
leading twist correction to the forward scattering amplitude of a gluon dipole Nq we substitute Eq. (93) 
into Eq. (52) obtaining 

Ng {r, b,y) = 2 Ni (r , b,y) + 2 N2 (r, b, y) - [N^ (r , hy)? + ■ ■ ■ , (94) 

where the first term on the right is the leading twist contribution Ng\ =2N\ given by Eq. (57), and the 
next two terms shown in Eq. (94) are the next-to-leading twist corrections. Higher twists are not shown in 
Eq. (94). To calculate the next-to- leading twist correction to gluon forward amplitude 

NG2{r,b,y) = 2N2{r,b,y) - [Ni{r,b,y)f = 2N2{r,b,y) - ^ [NGi{r,b,y)f (95) 

we use A^Gi from Eq. (57) and N2 calculated in [36]. Employing Eq. (23) from [36] in Eq. (9a) from the 
same reference would give us the first term on the right hand side of Eq. (95). In the end we write for a 
cylindrical nucleus 

J,^ r(^)r(^)r(i^A^) i . ,] 

\ r(i-^)r(i-f)r(Ai±Ai) 2[x(Ai) + x(A2)-x(Ai + A2)]^ 

The slight difference between the factors in the integrands of Eq. (96) and Eq. (23) of [36] is due to different 
definitions of the coefficients (cf. Eq. (15) of [36] with our Eq. (57)). 



^The argument presented in this paragraph is due to Larry McLerran. 
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FIG. 7. Gluon production in pA collisions at tlie next-to-leading twist level (see text). The blob in the center 
indicates a triple-pomeron vertex. 



The first term in the parentheses of Eq. (96) corresponds to the first term on the right hand side of 
Eq. (95). When we will substitute Ng2 from Eq. (96) into formula (56) for the cross section, this term 
would give the contribution illustrated in Fig. 7. It corresponds to the case when the gluon is produced still 
by the linear evolution with the triple pomeron vertex inserted below the emitted gluon. The rapidity of 
the triple pomeron vertex was integrated over in arriving at Eq. (96), with only the dominant contribution 
corresponding to the vertex being right next to the emitted gluon left [36] . (As was shown in [42] the diagrams 
where the triple pomeron vertex is inserted above the produced gluon cancel in the dipole evolution case 
considered here [37,35] in agreement with the expectation of the AGK cutting rules [66].) The second term 
in the parenthesis of Eq. (96) and on the right hand side of Eq. (95) corresponds to the case where the 
pomeron splitting occurs precisely at the rapidity position of the gluon production. The emitted gluon is 
produced by the first step of the non-linear evolution. (This term is the main difference between the results 
of [42] and [43].) As can be seen in the estimates performed below, this term contributes 50 ~ 100% of the 
answer depending on the /cT-region in question. 

Substituting Eq. (96) in Eq. (56) and integrating over z yields the following contribution to the gluon 
production cross section at the subleading twist level 



(Pk dy 



a, 2 (2^)3 J {2mf V kr 



CpSpSA f dAi c?A2 dA' ^ a „ ( Qso\ ^ / A 



+A2+A' _ Y 



g2a, a[x(Ai)+x(A2)]+2a, (y-y)x(A') 2Al+A2 + A'-3 V ^ J_r\ \\2 

Y (2 A1+A2+A' \ ^ 



A1+A2 



r(f)r(i 



Ai + As^ 
2 ^ 



r (1 - ^) r (1 - ^) r {^) 2[x{Xi) + x{M) - x(Ai + a^)] 



1 



(97) 



To study the onset of higher twist effects, we are interested in the next-to-leading twist contribution (97) 
in the region of transverse momenta /cT^fcgeom- Performing Ai and A2 integrations in Eq. (97) around the 
saddle point (64) and performing the A' integral in Eq. (97) around the saddle point (65) yields 



(Pk dy 



DLA 



as (27r)' 



lni/4 ^ 1^1/2 1^ 
(2a,y)3/2 (2S,(y-'y))3/4 




X exp 



4 W 2 tts y In 



QsO 



2d2as (Y-y) In 



A 



(98) 
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As one can see from Eq. (98) the next-to-leading twist correction tends to increase gluon production cross 
section at high-fcr- In the region of kr where Eq. (98) apphes, fcr^fcgeom, the higher twist corrections are 
paramctrically small and can not change the leading twist suppression of Eq. (74). 

To study higher twist effects in the extended geometric scaling region, Qs{y) ^ ^ ^geom, we evaluate Ai 
and A2 integrals in Eq. (97) around the saddle point (77) and do the A' integral around the saddle point (65) 
obtaining 



d'^k dy 



1^1/4 ^ 



LLA 



as (27r)5 



k^ 7C(3) a,y(2a,(y-y))3/4 



\2as{Y-y) 2 1nf^ 



In^ 



7C(3)a,y 



X exp 



2(ap-l)y + 2 J2a,(r-y) In^ ^^"^^ 



A 14C(3)a,t/ 



(99) 



The sign of Eq. (99) is dctcrniincxl by the sign of the expression in the parenthesis. One can see that for very 
large kr the expression in the parenthesis can become negative making the overall contribution to the cross 
section negative. However, Eq. (99) is valid only for kx < fcgeom and, therefore, can not be used at arbitrary 
high transverse momenta. At lower kr the sign changes and the term in Eq. (99) begins to contribute toward 
enhancement of RP^. The value of kx at which the sign transition takes place depends on the rapidity in 
question as well as on the atomic number A of the nucleus. To estimate the transition value of fcr one has 
to equate two terms in the parenthesis of Eq. (99). Assuming that In fcr/QsO ^ InQso/A we obtain 

A 

(100) 



with fcgeom given by Eq. (81). Therefore, the transition from suppression to enhancement in Eq. (99) happens 
at ko which is smaller than Ajgeom only by a factor of A~^/^^, which indicates that the term in Eq. (99) is pos- 
itive inside most of the extended geometric scaling region contributing to enhancement of gluon production. 
Here again one should note that Eq. (99) gives us a subleading twist contribution which is paramctrically 
smaller than the leading twist term from Eq. (86) in the kr region at hand (Qs(2/) ^ ^ fcgcom)- Eq. (99) 
is thus unlikely to affect the suppression of EF^ observed at the leading twist level in Eqs. (91) and (92). 

We conclude by observing that even after inclusion of quantum evolution (49) in the gluon production 
cross section, multiple rescatterings (higher twists) still tend to enhance gluon production at high-Zcy. In 
the kr region considered above, kr > Qs{y), these higher twist effects are still paramctrically small. In the 
next Subsection we are going to study the region of k-T where all twists become important, kx ~ Qs{y)- We 
will show that the combined effect of all twists is to introduce suppression of the Cronin maximum. 



D. Flattening of the Cronin peak 

We have demonstrated that the effect of quantum evolution (49) is to introduce suppression of RP^{kT, y) 
for kx^Qsiy) at the leading twist level. Let us now study what happens to RP"^{kT,y) at kx — Qs{y) as 
a result of evolution in energy. We showed in Sect. II that in the quasi-classical approximation the Cronin 
maximum of the ratio RP^(k, y) occurs at kx — QsO- In this Subsection we will follow the value of the ratio 
RP^{kT — Qs,y) to higher energies when quantum evolution is important. Since the position of the Cronin 
maximum is likely to be at kx — Qs{y) even when evolution is included, by studying R^^[kT = Qs^y) we 
are going to study the dependence of the height of the Cronin maximum on energy /rapidity. 

The fact that the scattering amplitude is a constant at the saturation scale [63,56,64] makes our calculation 
pretty straightforward. First we assume that Mellin transform of the gluon dipole amplitude obtained from 
the exact solution to the evolution equation Eq. (49) via Eq. (52) can be written as 



C^[zTQs{y)]\ ^T>^ 



geom 



(101) 
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The form of the solution presented in Eq. (101) assumes geometric scahng of No down to zt — l/Ajgeom 
and a leading twist expression for smaller zt in agreement with the analyses of [45,64]. Throughout this 
Subsection wc will use the definition of the saturation scale Qs{y) from Eq. (79) and the definition of /cgoom 
from Eq. (78). Our physical conclusions will be independent of the choice of definitions for saturation and 
geometric scales. 

Note that all information about the nonlinear evolution (49) is encoded in the fimction in Eq. (101). 
Using Eq. (101) in Eq. (56) we can calculate the differential gluon production cross section at kx = Qsiu)- 



Since fcgeom(y) ^ Qs(y) for large y we can set kg 



oo neglecting the Zt < 1/fcgeom part of the integral 



in Eq. (56). This approximation is justified in the Appendix A. We also assume that the dipole amplitude 
on a proton no is still given by the leading twist expression (60) around kx — Qsiv), which is a good 
approximation for a reasonable size nucleus. Substituting the first line of Eq. (101) and Eq. (60) into 
Eq. (56) we have 



cPkdy 



Cf Sp Sa A^ 



dX d\' 



dZTZTMQs{y)ZT) I ^^Cl,c^x^x'^ 



2'Ki 2Tri 



x(ztA)^'-2 [ztQM] 



A-2 „2a,{Y-y)x(y) 



(102) 



Performing the ^T-integration in Eq. (102) yields 
daP^ 



d?kdy 



kT=Qs{v) 



Cf Sp Sa f dX^ ^p ^a ^2 y2 ^ 
as7r(27r)2 } Im ^' ^ 



A+A'-3 



r(2-MV) 



— V e2«3(y-!/)x(A'). (103) 



Qs{y) 

It can be seen that all energy /rapidity and almost all atomic number dependence in Eq. (103) is given by the 
A'-integral. Since Qs{y) ^ A, the integral over A' in Eq. (103) can be evaluated in the double logarithmic 
approximation around the saddle point of Eq. (65) taken at /ct = Qs{y)- After that the integral over A 
carries almost no dynamical information. The result reads 



daP^ 



d^kdy 



where the integration over A gave an unknown function Ca defined as 



r(i-A)' 



Here we assume that Ca is only weakly (at most logarithmically) dependent on A, as is true for other 
coefficients like the one shown in Eq. (88). In case of an "ideal" geometric scaling the coefficient in 
Eq. (101) would be completely ^-independent ridding Ca of all of its A-dependence as well. 

To construct Rp^ we take the gluon production cross section in pp from Eq. (71) putting kr = Qs{y)- 
Substituting Eqs. (104) and (71) into Eq. (32) we get 



R'\Qs{y),y) = V^Ca {2a,f" ^ In^^ f ^^'l ^ exp ( -2 .^a.yln^) . (106) 



^+^T^ \ A J Q 



The energy and A-dcpcndence of Eq. (106) can be found using Eq. (79) and keeping in mind that Qso = 
A^/^ A. Since the definition of the saturation scale (79) is valid up to logarithmic prefactors, we can drop 
the prefactors in Eq. (106) leaving only 



BP^iQsivU) oc ^ exp (-2.j2asy\n^ ) . (107) 
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(If one defines Qs{y) by taking Nc{zt = ^/Qs{y),y) in the double logarithmic approximation and requiring 
that Ng{zt = 1/Qs{y), y) =const [56,63,64] the prefactors in Eq. (106) would cancel exactly.) Using Eq. (79) 
in Eq. (107) yields 



R''^{Qs{y),y) cx exp^4a,2/ ( l-^l + i^^ I } < 1. (108) 

We observe from Eq. (108) that in the course of quantum evolution the Cronin maximum of the ratio Rp^ 
decreases with energy until, at very high energy, it saturates at the lowest value RP"^ ~ A~^/^, which is 
much less than 1. The height of the Cronin peak is also a decreasing function of collision centrality/atomic 
number A, as can be seen from Eq. (108) 

Applicability of Eq. (108) is restricted by applicability of Eqs. (71) and (104). The latter two equations 
are valid only in the region where Qs{y) is larger than the geometric scale of the proton: Qs{y) > k^com/ A^^^ 
with fcgeom taken from Eq. (81), since this is where the transition between the saddle points takes place. 
With the help of Eq. (78) this condition becomes Qs{y) < A^^^QsO- In the kinematic region of extremely 
high y where this condition is not satisfied anymore, one has to replace Eqs. (104) and (71) by appropriate 
cross sections evaluated around the saddle point of Eq. (77). To generalize the conclusions presented above 
to arbitrary high rapidity let us follow Mueller and Triantafyllopoulos [64] and define the saturation scale 
by requiring that the power of the exponent in the leading twist expression for Nq given by (cf. Eq. (57)) 



J 27^^ 



is zero and stationary (its derivative with respect to A is zero) at zt = 1/Qs{y)- These conditions are satisfied 
at Ao = 1.255 [64] (our definition of x(A) is different from the one used in [64]). Resulting saturation scale is 
given by Eq. (83). Arguing that the A-intcgral in Eq. (57) is dominated by the saddle point in the cixponcnt we 
conclude that at zt = l/Qs{y) the gluon dipole amplitude Ng{zt = ^/Qs{y),y) is approximately constant 
at large y [64] 

f^C^ f TTTt) ^ ' = f^<^x e'"' ^'^^ HQMiy)) ^ constant(j/, A). (109) 

Making a similar assumption about the A'-integral in Eq. (103) taken at mid-rapidity (y = l"/2) and remem- 
bering that QsO = A^l^ K yields 



daP^ 



~ SaA-^"'^ ~ ^2/3-Ao/6_ (^^Q^ 

fcT=Qs(y),3/=y/2 



d^kdy 

Modifying Eq. (62) to give gluon production in pp at fcr = Qs{y) also taken at mid-rapidity we obtain 

..=Q.(.),.=./2 4a. (2.)3 j 2.^ 2.^ r (2 - A±v ) 



daPP 



d'^k dy 



A 



Qsiy) 



g2a3 [x(A)+x(A')]y_ 



Again the A- and A'-integrals in Eq. (Ill) are dominated by the saddle points at Aq giving an energy- 
independent cross section scaling as 



*We have recently learned that a similax conclusion regarding centrality dependence of the Cronin peak has been 
reached by A. Mueller and collaborators. 
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d?k dy 



~ (112) 

kT=Qs{y).y=Y/2 



with atomic number A. Combining Eqs. (110) and (112) with Eq. (32) yields 

R''^iQs{y),y) oc A-^/'^+^''/^ constant{y) ~ A'^-'^^^, (113) 

for Ao = 1.255. Note that the power of A in Eq. (113) is pretty close to that following from Eq. (108) and 
the two powers would be identical for Ao = 1. Note also that taking the expression for R^^ in the geometric 
scaling region from Eq. (92) and extrapolating it down to kr = Qs{y) one would obtain a power of A very 
close to that in Eq. (113) if one uses Qs{y) from Eq. (79). This conclusion not only verifies the self-consistency 
of our analysis, but also demonstrates that at asymptotic energies the height of Cronin maximum becomes 
(parametrically) equal to the height of the rest of the R^^ curve in the extended geometric scaling region. 
This is likely to indicate that at these energies the curve flattens out and the Cronin peak disappears. 

With the help of Eq. (113) we conclude that at high rapidities/energies the Cronin maximmn decreases 
with energy and centrality, with RP^{Qs{y),y) becoming less than 1. Eventually, at very high energy, the 
Cronin peak flattens out and saturates to an energy independent lower limit given by Eq. (113), which is 
parametrically suppressed by powers of A. 



E. Suppression Deep Inside Saturation Region 

Above we have shown that non-linear evolution (49) introduces suppression of gluon production in pA 
collisions making R^^ < 1 for kr^Qsiy)- In the region of smaller kr, kx <C Qs{y), we observed in Sect. 
IIIB that in the quasi-classical case of McLerran-Venugopalan model the ratio R^^ ^ 1 (see Eq. (47)). 
When the quantum evolution (49) is included it makes sense to consider the interval of low kx bounded 
from below by the saturation scale of the proton As{y), such that As{y) <C fcr *C Qs{y)- (For fcT< As(t/) 
the proton wave function also saturates and particle production in both pp and pA becomes similar to the 
case of AA, which has not been resolved even at the quasi-classical level [22,23,33]. Inclusion of evolution 
in A A is an even more difficult problem which we are not going to address here.) If kr is larger than the 
geometric scale of the proton fcgcomM^^^ (but still much less than Qsiy)) wc can use Eq. (71) to describe the 
gluon production cross section in pp. Deep inside the saturation region in pA the gluon production has been 
estimated in [42]. Employing Eq. (57) from [42] together with Eq. (71) we conclude that at mid-rapidity 



R'^iKeon. < kT « Qs{y),y) ~ ^ e-V2«..lnfc./A_ 



(114) 



sO 



Eq. (114) shows that inclusion of quantum evolution only introduces more suppression into R^^ at kr <C 

Qs iy) , making it a decreasing function of both the atomic number and energy. At very high energy k^ may 
become smaller than the geometric scale for the proton kgeom/A^^^ and the gluon production in pp would 
be driven by the saddle point (77) with A instead of Qso- Similarly to how it was done in [42] for DLA, 
one can estimate the gluon production cross section (56) deep inside the saturation region with the dipole 
amplitude on the proton evaluated around the LLA saddle point A' w 1. The result at mid-rapidity yields 



RP^iAsiy) « kr < kl 



geomy 



^-1/3 



Qs{y) 
A 



exp 



-{ap - l)y + 



2\n\kTlK)~\n\Q,{y)/h) 



14C(3)a.y 



(115) 



Therefore, at very high energies the ratio R^^ becomes almost independent of kx even at very low kx- Using 
the saturation scale from Eq. (79) in Eq. (115) at asymptotic energies gives 



RP^{K{y)<^kx<kP) ~ A-o-^e 



-1,0 OLs y 



(116) 



We observe again that nonlinear evolution leaves Rp^ very small at kx < QAjj)- R^^ given by Eq. (116) 
is a decreasing function of both rapidity/energy and centrality. This conclusion seems natural, since the 
saturation effects are known to soften the low-fcr gluon spectra in pA compared to pp. 
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F. Toy Model 



To illustrate the conclusions reached above let us construct a simple toy model exhibiting suppression of 
R^^ at all kr- We start with the quasi-classical formula for gluon production in pA in the following form 
which could be obtained from Eq. (40) for a cylindrical nucleus and for azimuthally symmetric Nq 



(Pk dy 



cisCf Sa 
7r2 fc2 



f 

Jo 



1 



dzT Joikr zt) In ^-^ [zt d;,^NG{zT, y = 0)]. 



(117) 



With the increase of energy the gluon dipolc amplitude on the nucleus will reach saturation. Therefore, its 
2;T-dependence will change more significantly than for the corresponding amplitude on the proton, which 
will stay unsaturated. (Of course at very high energy the dipole amplitude on the proton will also reach 
saturation, but we are not going to consider that energy range here.) Therefore, in our toy model we will 
assume for simplicity that the gluon dipole amplitude on the proton remains unchanged with increasing 
energy, giving Inl/ [zt^^) in Eq. (117). We will model the gluon dipole amplitude at high energy by a 
Glauber-like unitary expression 



1 



-^TQs{y) 



(118) 



which mimics the onset of anomalous dimension A = 1 by the linear term in the exponent. The saturation 
scale Qs{y) in Eq. (118) is some increasing function of y which can be taken from Eq. (79) or from Eq. (83). 

Indeed the amplitude in Eq. (118) has an incorrect small-zy behavior, scaling proportionally to zt instead 
of z"^ as shown in Eq. (14). If Eq. (118) is used in Eq. (117) it would lead to an incorrect high-fcr behavior 
of the resulting cross section. We therefore argue that Eq. (118) is, probably, a reasonable model for Nq 
inside the saturation and extended geometric scaling regions {I/zt ~ /ct < Aigeom), but should not be used 
for very small Zt / high kr {1 / Zt kr > fcgeom)- 

Substituting Eq. (118) into Eq. (117) and integrating over zt yields 



da\ 



pA 
toy 



d?k dy 



krp krp -\- 



Qs{kl + Ql) + ^kl + Ql 



In 



2(fc|, 



krA 



+^ In 



Qs 



2 ^/W+Ql + Qs 



(119) 



where 7 is the Euler's constant and Qs = Qs{y)- Corresponding gluon production cross section for pp is 
obtained by expanding Eq. (119) to the lowest order at high kr and substituting A instead of Qs and Sp 
instead of Sa- 



d?k dy 



cusCf SpA 



In — h 7 



(120) 



Of course in Eq. (120) one implicitly assumes that anomalous dimension has set in for only one of the 
protons in pp. This assumption is not valid at mid-rapidity, but may be used to study particle production 
at rapidities near fragmentation region of one of the protons. 
Substituting Eqs. (119) and (120) in Eq. (32) yields 



moy{kT,y) 



krA 



Qs (fc2 +g2)(in(2fcy/A)+7) 



■ <3«(4 + Qs) + \Jkl + Ql {^Ql + 7 



1^2 2(fcj + Ql) ^ "s/^T + ~ 
^ krA 2 ^fc^ + Q1 + Qs 



(121) 



in which we assumed that A is the saturation scale of the proton such that Ql = A^/^ A^ even at high energy. 
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FIG. 8. The ratio i?*"* plotted as a function of kr/Qs for (i) McLcrran-Venugopalan model, which is valid for 
moderate energies (upper solid line); (ii) our toy model for very high energies/rapidities from Eq. (121) (lower solid 
line); (iii) an interpolation to intermediate energies (dash-dotted and dashed lines). The cutoff is A = 0.3 Qs- 

The toy model i?f^(A:T,y) from Eq. (121) is plotted as a function of kr/Qs in Fig- 8 for A = 0.3 Qs 
(lower solid curve). It exhibits suppression of gluon production in pA at all values of kr leveling off at 
R^^y ^ A/Qs ~ A~^l^ for kr^Qs at high energy, in agreement with our conclusions of Sections IIIB and 
HID. 

Our toy model (121) represents the high energy asymptotics of . To compare it to lower energies, we 
also plot BP^ for the quasi-classical McLerran-Venugopalan model given by Eq. (34) (upper solid curve in 
Fig. 8). As the energy increases the upper solid line in Fig. 8 would decrease eventually turning into the 
lower solid line. The corresponding intermediate energy stages are shown by the dash-dotted and dashed 
lines in Fig. 8. These lines are for illustrative purposes only and do not correspond to any toy model. They 
demonstrate how the Cronin peak gradually disappears as energy or rapidity increase. 



V. CONCLUSIONS 



In this paper we have demonstrated that saturation effects in the gluon production in pA at moderate 
energy can be taken into account in the quasi-classical framework of McLerran-Venugopalan model, which 
includes Glauber-Mueller multiple rescatterings, resulting only in Cronin enhancement of produced gluons 
at fcr = (1 2) Qs0 7 as was shown in Fig. 4 and in Eq. (37). Similar conclusions have been reached in [48]. 
In this quasi-classical approximation the height and position of the Cronin peak are increasing functions of 
centrality as indicated by Eq. (38). 

We have also shown that at higher energies/rapidities, when quantum evolution becomes important, it 
introduces suppression of gluons produced in pA collisions at all values of fcr, as compared to the number 
of gluons produced in pp collisions scaled up by the number of collisions Ncoii, as suggested previously [8]. 
The resulting RP^ at high energy /rapidity is a decreasing function of centrality. We have considered three 
different complimentary regions of kr, which cover together all of fcr-range: 

i. kT > Qs{y) region. Gluon production cross section in pA is dominated by the leading twist effects in 
this region of kr- We have shown how the leading twist suppression arises in the double logarithmic 
approximation for kr > fcgeom '> Qsiv) with the corresponding RP-^{kT,y) given by Eq. (72), which 
approaches 1 as A:t ^ oo- At Qs{y) < ^T^^goom the leading twist suppression is mainly due to the 
change in anomalous dimension A from its double logarithmic value (64) to the leading logarithmic 
value (77) . R^^ {kx , y) for this kx-'wmdow is given by Eq. (89) leading to suppression described by 
Eq. (91). At very high energies, when the extended geometric scaling regions of the proton and the 
nucleus overlap (for Qs{y) < kr^ ^|eom) the decrease of R^"^ with energy stops at roughly R^"^ ^ A~^/^ 
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as follows from Eq. (92). This leading twist effect has been originally pointed out in [8]. We have not 
considered suppression mechanisms that may stem from running of the coupling constant, which would 
modify the A-dependence of the saturation scale [65]. 

ii. kx ^ Qs{y) is the position of the Cronin maxinmm in the quasi-classical approximation. We began 
the analysis of this fcr-region by studying higher twists in the adjacent region of kr > Qs{y)- Next- 
to-leading twist term was shown to contribute towards enhancement of R^"^ at high-fcr even when 
evolution is included. However, higher twist effects are parametrically small at fc^ > Qs{y) and can 
not change our leading twist conclusions about suppression. To assess the contribution of all twists 
we studied the behavior of the Cronin maximum (/ct ~ Qs{y}) with increasing energy. We showed 
that at kr = Qs{y) is a decreasing function of energy /rapidity and centrality saturating at the 
energy-independent lower bound given by Eq. (113). Since the height of the Cronin maximum becomes 
parametrically of the same order as the rest of RP"^ at higher kr given by Eq. (92), we conclude that 
Cronin peak disappears at asymptotically high energies/rapidities. 

iii. kx "C Qs {y) region. The suppression of R^^ deep inside the saturation region, kx "C Qs {y) , only gets 
stronger as the evolution (49) is included (see Eq. (116)). 

Our results are summarized in Fig. 8. 

It is interesting to observe that the behavior of R^^ at high energies is qualitatively different from what 
one would expect by taking the quasi-classical expression (34) and letting Qs in it increase with energy. In 
case of DIS a similar trick where one replaces Qso in the Glauber- Mueller expression for the dipole cross 
section (9) by the energy dependent Qs from, for instance, Eq. (79) leads to correct qualitative behavior 
of resulting F2 structure function and even generates some successful phenomenology [67]. However, as we 
showed above, a naive generalization of McLerran-Venugopalan model by increasing Qg with energy does 
not work for Rp^ even at the qualitative level. 

The analysis in the paper was, of course, done for sufficiently high energy and/or rapidity, such that the 
saturation approach was assumed to be still valid for the highest kr involved. This implies that the effective 
Bjorken x is still sufficiently small for all kx we consider. The extent to which this treatment applies at 
high kr hadron production at RHIC is difficult to assess theoretically. We thus eagerly await the results of 
the experimental analyses of centrality dependence of hadron production above the Cronin region (kx > 6 
GeV). It is also very important to extend the present measurements away from the central rapidity region 
to separate initial state effects from possible energy loss in cold nuclear matter. Indeed, in the deuteron 
fragmentation region, the effects of saturation in the Au wave function will be enhanced, while the density of 
the produced particles (see, e.g., the predictions in [68]) and thus the associated energy loss will be minimal. 
In the A21 fragmentation region the opposite will be true. 

We, therefore, conclude that if the effects of quantum evolution and anomalous dimension are observed in 
the forward rapidity region of dAu collisions at RHIC, they would manifest themselves by reducing R'''^ at 
all kx as shown in Fig. 8, eliminating the Cronin enhancement. R'^^ will become a decreasing function of 
centrality. The pA program at LHC would observe an even stronger suppression of Rp^. However, it might 
be that the quantum evolution effects are still not important even in the forward region of dAu collisions at 
RHIC. Then reduction of R'^^ going from mid-rapidity to deuteron fragmentation region should be rather 
mild and the Cronin peak would not disappear in the forward region. The relevant particle production 
physics would be described by McLerran-Venugopalan model. The height of the Cronin peak would then be 
an increasing function of centrality. 

If the forthcoming data on i?*^"* in the forward rapidity region of dAu collisions would have no high-pT 
suppression and would exhibit only a strong Cronin maximum which is an increasing function of centrality 
in agreement with predictions of multiple rescattering models described in Sect. Ill [48,49,51-54], then all of 
the observed high-pr suppression in Au — Au collisions would have to be attributed to the final state effects. 
However, if the future R"^^ data in the forward rapidity region exhibits suppression either for all px or at 
high Px with R"^^ being a decreasing function of centrality as described in this paper (see also [8]), then 
a fraction of R^^ suppression in the forward rapidity region of Au — Au collisions should be attributed to 
initial state quantum evolution effects. Indeed, there is some evidence [4] that the high kx suppression in 
Au — Au collisions increases between the pseudo-rapidities 77 = and rj = 2.2. 

The dAu data at y ~ [1-4] also suggest suppression of the yields of charged hadrons [3] and neutral 
pious [1] at kx > 6 GeV, though the suppression is not significant statistically. If this initial-state effect is 
confirmed, it should also be taken into account in the interpretation of Au — Au results at ?/ ~ 0. 



28 



The dAu results will thus allow to clarify the relative importance of initial and final state interactions 
at different transverse momenta and rapidities of the produced particles. They will be indispensable for 
establishing a complete physical picture of heavy ion collisions at RHIC energies. 

Note added: After the first version of this paper appeared, a similar analysis has been done in [69-71]. 
The analyses of [69-71] agree with our conclusions on the presence of Cronin effect in the quasi-classical 
approximation. The results of [69,71] are also in agreement with our conclusion about high-pr suppression 
of gluon production. 
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APPENDIX A 



Here we are going to derive Eq. (102). In writing down Eq. (102) we assumed that the integration region 
Zt < l/fcgeom is negligible. To justify this approximation let us start by substituting Eq. (101) into Eq. (56). 
We find 



daP^ 



d^kdy 



dZT ZT Jo{Qs{y)zT) 



kT=Qs(y) 



CpSpSA [ 

a,^(2^)2g2(y) 
X (^T A)^'-2 (^T g.0)^"' y x(A)+2 a, (Y -y) x(A') 



27ri ./ 2m ^ ^ ^ 



+ 



j^^ dzT ZT MQs{y)zT) I ^J'^iCl, X' X'^ A' Qliy) 



1/fcgeom 

x{zTAf-' [zTQs{yt-'e'^'^''-'^>^^^'^y 
The difference between Eq. (Al) and the target Eq. (102) is 



Cf Sp Sa 
a,7r(27r)2Q2(y) 



I'^'-dzrz^^JomvW) I ^J'^,Cl,X^X'^ 



sO 



(Al) 



{ZT Qso)^ e2«^ ^>^(^) - [zT Qs{y)]^ 



(A2) 



Since fcgeom ^ Qs{y) we can neglect the argument of the Bessel function in the integral in Eq. (A2) putting 
•^o(O) = 1. Integration over zt then yields 



Cf Sp Sa f dX f dX' p .2 „2as{Y-y)x{y) 
a.7r(27r)2 J 2ni J 2ni ^ ^ ^ 



A 



A' 



A + A' 



'"geom 



Q 



sO 



(A3) 
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Due to the inequality fcgeom ^ Qs{y) ^ A, the integration over A' in Eq. (A3) is dominated by the saddle 
point at A' « 2, as shown in Eqs. (63) and (65). The integral over A in Eq. (A3) becomes 



/ 

f dX 
J 2^ 



dX X^ 
X + X' -2 

dX 



\ ^geom / 



A 



Qs{y) 



X 



- C< 



h I ^ \ h 

"'geom / \ "'geom 



d 



9 In k. 



gcom 



[Ng{zt ^ (l/fcgeom)",y) " Ng{zt ^ (lAgeom)^ , y)] = 0, 



(A4) 



where we assumed that Ng{z, y) and its derivatives with respect to zt from Eq. (101) are smooth functions 
of Zt such that the difference of the above limits is zero. This assumption is justified since NG{z,y) is 
proportional to the scattering matrix which is an analytic function of its variables. Eq. (49) makes Ng{z, y) 
analytic by construction. 

We showed that the difference between the exact Eq. (Al) and our Eq. (102) is zero, making the two 
equations equal, as desired. However, the above proof required that the representation of NG{z,y) given by 
Eq. (101) has a smooth matching of the two regions at Zt = 1/fcgeom, i-e., that representation (101) is not 
just a good approximation but an exact identity. To show that no such assumption is required to prove that 
the expression in Eq. (A3) is a negligible correction to Eq. (102) let us estimate the energy dependence the 
first term in Eq. (A3). The second term in Eq. (A3) is negative and can only make the overall contribution 
smaller. Employing double logarithmic approximation for A- and A'-integrals and using Eqs. (79) and (78) 
we derive (setting y = Y/2 for simplicity) 



dX 



dX' 



Cl, A^A'^ e2«3(y-?/)x(V) 



A + A' 



A 



A' 



A-2 

"gcom 



,2-cisyxW 



Q 



sO 



^seom Qsiy) 



DC A 



'6V2 e 



-4as 2/(3-2^2) 



(A5) 



which is a decreasing function of rapidity and ccntrality. It is obviously negligible compared to the increasing 
function of y given by Eq. (104). This accomplishes our proof of Eq. (102). 
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